{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "cell_id": "00000-b62acb5e-bf38-4990-b0b2-de1a4ccf4944",
    "output_cleared": false
   },
   "source": [
    "# SPECT simulation \n",
    "\n",
    "Analysis after the simulation :\n",
    "- in folder: spect/\n",
    "- macro: main.mac"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {
    "cell_id": "00001-ec633223-4fbf-496a-990e-5c5826265081",
    "execution_millis": 2,
    "execution_start": 1605009142403,
    "output_cleared": false,
    "source_hash": "206465ed"
   },
   "outputs": [],
   "source": [
    "%matplotlib widget\n",
    "import matplotlib.pyplot as plt\n",
    "from matplotlib import collections  as mc\n",
    "import numpy as np\n",
    "from ipywidgets import interact\n",
    "import os\n",
    "from pathlib import Path\n",
    "import uproot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "cell_id": "00002-5be4f8b1-56ff-4026-8ab0-0934df8e2854",
    "execution_millis": 35,
    "execution_start": 1605009142412,
    "output_cleared": false,
    "source_hash": "bd411e75"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The Current Working Directory (CWD) is: \n",
      " /Users/dsarrut/src/gate/dqprm/2021/gate-exercices/spect\n"
     ]
    }
   ],
   "source": [
    "# The following command display the current working directory (where jupyter has been launched)\n",
    "cwd = os.getcwd()\n",
    "print('The Current Working Directory (CWD) is: \\n', cwd)\n",
    "folder = Path('./output_ref')                            # <--------------------------------------TO CHANGE BY YOUR OWN OUTPUT FOLDER \n",
    "if (not folder.is_dir()):\n",
    "    print('ERROR: {} is not a folder.'.format(folder))\n",
    "filename = os.path.join(folder, 'spect.root')\n",
    "f = uproot.open(filename)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "cell_id": "00003-c0a724eb-32bf-4cef-94c4-9a9d5913e516",
    "execution_millis": 4,
    "execution_start": 1605009142450,
    "output_cleared": false,
    "source_hash": "365ca338"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "List of keys: \n",
      " [b'latest_event_ID;1', b'total_nb_primaries;1', b'Hits;1', b'SinglesAdder;1', b'SinglesBlurring;1', b'SinglesThresholder;1', b'SinglesUpholder;1', b'Singles;1', b'peak171;1', b'scatter;1', b'peak245;1']\n"
     ]
    }
   ],
   "source": [
    "print(\"List of keys: \\n\", f.keys())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "cell_id": "00004-740bc229-33ef-4fb7-ae32-98d6b3c28afb",
    "execution_millis": 61,
    "execution_start": 1605009142457,
    "output_cleared": false,
    "source_hash": "a77506f"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of hits: 22007\n",
      "Hits: \n",
      "b'PDGEncoding'  b'trackID'  b'parentID'  b'trackLocalTime'  b'time'  b'edep'  b'stepLength'  b'trackLength'  b'posX'  b'posY'  b'posZ'  b'localPosX'  b'localPosY'  b'localPosZ'  b'momDirX'  b'momDirY'  b'momDirZ'  b'headID'  b'crystalID'  b'pixelID'  b'unused3ID'  b'unused4ID'  b'unused5ID'  b'photonID'  b'nPhantomCompton'  b'nCrystalCompton'  b'nPhantomRayleigh'  b'nCrystalRayleigh'  b'primaryID'  b'sourcePosX'  b'sourcePosY'  b'sourcePosZ'  b'sourceID'  b'eventID'  b'runID'  b'axialPos'  b'rotationAngle'  b'volumeID'  b'processName'  b'comptVolName'  b'RayleighVolName'  "
     ]
    }
   ],
   "source": [
    "hits = f[b'Hits']\n",
    "print('Number of hits:', len(hits))\n",
    "print('Hits: ')\n",
    "for i in hits.keys():\n",
    "    print(i, ' ', end='')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "cell_id": "00005-8b47ee00-7104-4efa-9801-0b2fe98e45d5",
    "output_cleared": false
   },
   "source": [
    "# Plot 1: position of hits, hits times (global & local)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {
    "cell_id": "00006-a1153054-743a-4484-b3d7-256edea767c7",
    "execution_millis": 3252,
    "execution_start": 1605009142519,
    "output_cleared": false,
    "source_hash": "7f8e5400"
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "bca14deefd6f4e148ddbde99f9a421b1",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "Text(0.5, 0, 'Time in ns')"
      ]
     },
     "execution_count": 39,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(1, 3, figsize=(12, 4))\n",
    "\n",
    "# Plot 1\n",
    "hits_u = hits.arrays()[b'posX']\n",
    "hits_v = hits.arrays()[b'posZ']\n",
    "ax[0].scatter(hits_u, hits_v, alpha=0.8, s=3)\n",
    "ax[0].set_title('Positions (2D) of hits in the camera')\n",
    "ax[0].set_xlabel('X position in mm')\n",
    "ax[0].set_ylabel('Y position in mm')\n",
    "#ax[0].set_aspect('equal')\n",
    "\n",
    "# Plot 2\n",
    "times = hits.arrays()[b'time']\n",
    "n, bins, patches = ax[1].hist(times, 100, facecolor='green', alpha=0.75)\n",
    "ax[1].set_title('Hits time (global)')\n",
    "ax[1].set_xlabel('Time in s')\n",
    "\n",
    "# Plot 3\n",
    "times = hits.arrays()[b'trackLocalTime']*1e9\n",
    "n, bins, patches = ax[2].hist(times, 100, facecolor='green', alpha=0.75)\n",
    "ax[2].set_title('Local hit times (time from source emission)')\n",
    "ax[2].set_xlabel('Time in ns')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "cell_id": "00007-d230344f-c2fe-4569-be4d-1a75fdad9663",
    "output_cleared": false
   },
   "source": [
    "# Plot 2: edep, step length and track length"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {
    "cell_id": "00008-e2b8bf16-0053-4df6-85a9-00bdca67e379",
    "execution_millis": 2701,
    "execution_start": 1605009145833,
    "output_cleared": false,
    "source_hash": "fb76ead1"
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "767833c16d1a4a41b0383799f42420a2",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "Text(0.5, 0, 'Length in mm')"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fig2, ax2 = plt.subplots(1, 3, figsize=(12, 4))\n",
    "\n",
    "# Plot 1\n",
    "edep = hits.arrays()[b'edep']\n",
    "n, bins, patches = ax2[0].hist(edep, 100, facecolor='green', alpha=0.75)\n",
    "ax2[0].set_title('Edep')\n",
    "ax2[0].set_xlabel('Energy in MeV')\n",
    "\n",
    "# Plot 1\n",
    "stepLength = hits.arrays()[b'stepLength']\n",
    "n, bins, patches = ax2[1].hist(stepLength, 100, facecolor='green', alpha=0.75)\n",
    "ax2[1].set_title('Step Length')\n",
    "ax2[1].set_xlabel('Length in mm')\n",
    "\n",
    "# Plot 2\n",
    "trackLength = hits.arrays()[b'trackLength']\n",
    "n, bins, patches = ax2[2].hist(trackLength, 100, facecolor='green', alpha=0.75)\n",
    "ax2[2].set_title('Track Length')\n",
    "ax2[2].set_xlabel('Length in mm')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {
    "cell_id": "00009-e878ae6a-d878-4821-8fde-7698918ba06a",
    "execution_millis": 7,
    "execution_start": 1605009148550,
    "output_cleared": false,
    "source_hash": "10778042"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of Singles : 14545\n",
      "Singles: \n",
      "b'runID'  b'eventID'  b'sourceID'  b'sourcePosX'  b'sourcePosY'  b'sourcePosZ'  b'time'  b'energy'  b'globalPosX'  b'globalPosY'  b'globalPosZ'  b'headID'  b'crystalID'  b'pixelID'  b'unused3ID'  b'unused4ID'  b'unused5ID'  b'comptonPhantom'  b'comptonCrystal'  b'RayleighPhantom'  b'RayleighCrystal'  b'axialPos'  b'rotationAngle'  b'comptVolName'  b'RayleighVolName'  "
     ]
    }
   ],
   "source": [
    "singles = f[b'Singles']\n",
    "print('Number of Singles :', len(singles))\n",
    "print('Singles: ')\n",
    "for i in singles.keys():\n",
    "    print(i, ' ', end='')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "cell_id": "00010-6b9ee439-23f5-4451-be31-751890e44398",
    "output_cleared": false
   },
   "source": [
    "# Plot 3: energy  of the singles"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {
    "cell_id": "00011-e7c5fc9c-0e27-4bb6-8821-438cf7653ff5",
    "execution_millis": 577,
    "execution_start": 1605009148571,
    "output_cleared": false,
    "source_hash": "401d2dd6"
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "9d558d7ae3234628a07550b9217e5a4f",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "Text(0.5, 0, 'Energy in keV')"
      ]
     },
     "execution_count": 42,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fig3, ax3 = plt.subplots()\n",
    "energy = singles.arrays()[b'energy']/0.001\n",
    "n, bins, patches = ax3.hist(energy, 200, facecolor='green', alpha=0.75)\n",
    "ax3.set_title('Edep')\n",
    "ax3.set_xlabel('Energy in keV')"
   ]
  }
 ],
 "metadata": {
  "deepnote_execution_queue": [],
  "deepnote_notebook_id": "2fdaf907-cc5d-4f6c-96dc-43e4d6245abf",
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
